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1.  Introduction 


This  report  defines  algorithms  that  can  be  used  to  calculate  speed,  distance,  and  time  variables  as 
functions  of  initial  and/or  final  projectile  states.  Examples  include  speed  as  a  function  of  initial 
speed  and  distance  traveled,  distance  traveled  as  a  function  of  initial  and  final  speeds,  and  time  of 
flight  as  a  function  of  initial  and  final  speeds.  Derivations  for  the  algorithms  are  based 
exclusively  on  force  due  to  air  resistance  and  allow  for  drag  coefficients  that  have  a  functional 
relationship  to  Mach  number. 

All  algorithms  are  written  in  terms  of  scaled  variables  (speed,  distance,  and  time)  that  are 
independent  of  units.  This  approach  allows  for  the  precalculation  of  a  set  of  parameters  that 
depend  only  on  the  characteristics  of  a  particular  drag  versus  Mach  curve  (i.e.,  the  parameters  are 
independent  of  projectile  mass,  projectile  cross-sectional  area,  and  air  density).  A  practical 
benefit  of  this  approach  is  a  reduction  in  computation  time. 

The  algorithms  presented  in  this  report  are  particularly  useful  for  situations  where  total  distance 
traveled  is  approximately  equal  to  the  distance  between  starting  and  ending  locations.  Examples 
include  flat-fire  trajectories  and  fragment  trajectories. 

JTCG/ME-79-1-2*  presents  a  method  for  calculating  speed  as  a  function  of  initial  speed  and  total 
distance  traveled.  It  relies  upon  an  analytic  solution  that  is  similar  to  what  is  presented  in 
sections  3.1  and  3.2  of  this  report.  The  key  difference  is  that  the  method  presented  in  this  report 
provides  a  solution  to  the  equations  of  motion  in  terms  of  variables  that  are  independent  of  units. 


2.  Scaled  Variables 


Recall  that  Mach  number,  M,  is  defined  to  be  the  ratio  of  the  speed,  v,  of  an  object  moving 
through  a  fluid  to  the  local  speed  of  sound,  v^.  That  is, 

(1) 

Scaled  distance,  is  defined  as  a  function  of  measured  distance,  x: 


'joint  Technical  Coordinating  Group  for  Munitions  Effectiveness.  Computer  Program  for  General  Full  Spray  Materiel  MAE 
Computations,  Volume  II  -  Analyst  Manual;  AMSAA  61  JTCG/ME-79-1-2;  U.S.  Army  Materiel  Systems  Analysis  Aetivity: 
Aberdeen,  MD,  1987. 
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(2) 


2m 


where  p  is  the  density  of  the  air  through  which  a  projectile  travels,  A  is  the  projectile’s  cross- 
sectional  area,  and  m  is  its  mass. 


Scaled  time,  z ,  is  defined  as  a  function  of  measured  time,  t: 


T  = 


(3) 


3.  Derivation  of  Equations  for  Linear  Drag  Coefficients 


Begin  with  the  equation  for  the  force  acting  on  a  projectile  that  is  due  to  air  resistance: 

F  =  -^C^pAv^  (4) 

where  Co  is  a  projectile’s  drag  coefficient,  and  p ,  A,  and  v  are  as  defined  in  section  2. 

3.1  Drag  Coefficient  as  a  Function  of  Mach  Number 

Assume  that  Cd  is  linearly  dependent  on  Mach  number  (figure  1).  Then 


C^=aM  +  l3  (5) 

where  a  and  p  are  constants  representing  slope  and  intercept,  respectively. 
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3.2  Scaled  Distance  as  a  Function  of  Mach  Number 

Using  Newton’s  second  law,  equation  4  can  be  used  to  solve  for  acceleration: 

F  =  ma  (6) 

^  a  =  -^Cj^pAv^  (7) 

2m 


Next,  the  definitions  of  speed  and  acceleration  can  be  used  to  find  acceleration  as  a  function  of 
speed  and  position: 


_  dx 
dt 


,  dv  dv  dx 

and  a  =  —  = - 

dt  dx  dt 


dv 

=>  a  =  V  — 
dx 


(8) 

(9) 

(10) 


Combining  equations  7  and  10, 


dv 

dx 


-^CdPAv 

2m 


(11) 


Next,  equations  1  and  2  can  be  used  to  rewrite  equation  1 1  in  terms  of  scaled  variables: 


dM 

ds 


=  -CoM 


(12) 


Substituting  equation  5  into  equation  12, 

dM 

ds 


=  -{aM  +  P)M 


Equation  13  can  be  solved  for  s  by  using  the  method  of  separation  of  variables: 

1 


=  -l; 


■dM 


{aM  +  P)M 

The  solution  to  the  integral  in  equation  14  can  be  found  in  integral  tables:^ 


s  =  —In 

P 


ct  + 


l_ 

M 


+  C 


(13) 


(14) 


(15) 


where  c  is  a  constant  of  the  integration  and  is  determined  by  initial  values  for  s  and  M. 


2 

Lide,  D.  R.  Handbook  of  Chemistry  and  Physics,  89th  ed.;  CRC  Press:  London,  2008;  p  A-16,  equation  37. 
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(16) 


c-  Sr, - In 

P 


oc  + 


p 


M, 


0/ 


Note  that  equations  15  and  16  are  only  valid  for  a  +  pi  M  >(),  M  ^0,  and  p  ^0.  Since,  by 
definition,  Mach  numbers  are  non-negative,  the  conditions  a  +  p/M  >0  and  M  ^0  are  assured 
if  we  assume  that  >  0  for  all  values  of  M  and  >  0.  The  case  where  p  =  0  will  occur 
occasionally.  To  be  exact,  equation  14  should  be  solved  for  the  special  case  where  p  =  0  An 
practice,  however,  when  encountering  a  situation  where  =  0 ,  it  is  sufficient  to  change  p  to 
some  very  small  value. 

Figure  2  presents  scaled  distance,  s,  as  a  function  of  Mach  number,  M. 


Figure  2.  Scaled  distance,  s,  as  a  function  of  Mach  number,  M. 


3.3  Mach  Number  as  a  Function  of  Scaled  Distance 

Solving  equation  15  for  M  yields  Mach  number  as  a  function  of  scaled  distance  (figure  3): 


(17) 
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3.4  Scaled  Time  as  a  Function  of  Mach  Number 


Begin  by  using  equations  1  and  7  to  find  acceleration  as  a  function  of  Mach  number: 


a  =  -^C^pAvlM^ 

2m 

Next,  use  equations  1  and  9  to  relate  acceleration  to  the  time  derivative  of  Mach  number: 

dM 


a  =  v. 


dt 


Combining  equations  18  and  19, 


dM  1 


C^pAvX^ 


dt  2m 

Equation  3  can  be  used  to  rewrite  equation  20  in  terms  of  scaled  time: 

dM 


dr 


(18) 


(19) 


(20) 


(21) 


Substituting  equation  5  into  equation  2 1 , 

dM 

dr 


=  -(aM  +  J3)M^ 


Equation  22  can  be  solved  for  r  by  using  the  method  of  separation  of  variables: 

1 


■dM 


{aM  +  P)M^ 

The  solution  to  the  integral  in  equation  23  can  be  found  in  integral  tables:^ 


a 


T  =■ 


1 


/3M 


In 


cc  +  ■ 


M 


+  d 


(22) 


(23) 


(24) 


where  d  is  a  constant  of  the  integration  and  is  determined  by  initial  values  for  r  and  M  . 


d  =  Tq 


1 

Wl 


a  .  , 


( 


OC  + 


V 


0  7 


(25) 


As  was  the  case  with  equations  15  and  16,  equations  24  and  25  are  only  valid  for  a  +  p/M  >0, 
M  and  p^O.  Since,  by  definition,  Mach  numbers  are  non-negative,  the  conditions 
a  +  pi  M  >  0  and  M  ^0  are  assured  if  we  assume  that  >  0  for  all  values  of  M  and  >  0. 


Lide,  D.  R.  Handbook  of  Chemistry  and  Physics,  89th  ed.;  CRC  Press:  London,  2008;  pp  A- 17,  equation  40. 
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The  case  where  /3  =  0  will  occur  occasionally.  To  be  exact,  equation  23  should  be  solved  for  the 
special  case  where  =  0 .  In  practice,  however,  when  encountering  a  situation  where  =  0 ,  it  is 
sufficient  to  change  p  to  some  very  small  value. 


Figure  4  presents  scaled  time,  r  ,  as  a  function  of  Mach  number,  M. 


4.  Generalization  to  an  Arbitrary  Drag-Mach  Relationship 


The  equations  in  section  3  can  be  generalized  for  a  piecewise,  linear  function  of  drag  given 
Mach.  This  technique  allows  for  any  functional  relationship  between  drag  coefficient  and  Mach 
number  to  be  modeled. 

4.1  Drag  Coefficient  as  a  Piecewise,  Linear  Function  of  Mach  Number 

Suppose  that  drag  coefficient  as  a  function  of  Mach  number  is  defined  by  n  matched  pairs  of 
values  for  Mach  number  and  drag  coefficient,  (M^ ,  ) ,  such  that 

for  0<k<n-\  (26) 

Using  the  matched  pairs,  n  - 1  equations  that  have  the  form  of  equation  5  can  be  generated: 

C^=a,M  +  fi,  (27) 


for  k  such  that  M  <M<  and  0<k  <n-\. 

Figure  5  presents  drag  coefficient,  Cd,  as  a  piecewise  function  of  Mach  number,  M. 
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Equation  27  can  be  used  to  determine  the  values  of  and  : 


C  -  C 

^D,k  '^D,k+l 

ar  = - 

and 

~  ^ D,k  ~  ^k^ k 

4.2  Scaled  Distance  as  a  Function  of  Mach  Number 

Equation  15  can  be  generalized  to  a  piecewise  function  of  Mach  number: 

1 


5  = - In 

Pk 


a. 


A 

M 


+  Cl 


for  k  such  that  <M<  M . 


Values  for  can  be  found  by  generalizing  equation  16: 


^k  ~  ^Q,k  r. 

Pk 


Uk  + 


Pk 


M, 


(28) 


(29) 


(30) 


(31) 


V  o,i  J 

Values  for  initial  conditions  for  each  segment  of  the  function  are  given  by  equations  32  and  33: 


M,,,  =M, 


(32) 


and 

So,k^s,  (33) 

Values  for  s^.  can  be  obtained  from  equation  30: 


7 


1 


(34) 


A-i 


In 


a 


A-i 


k-\ 


M 


+  c 


k-i 


k  J 


Since  the  value  of  is  arbitrary,  zero  is  used  as  a  convenient  choice: 


^^=0  =  0 


(35) 


Figure  6  presents  an  overview  of  the  functions  necessary  to  calculate  scaled  distance  as  a 
function  of  Mach  number.  Note  that  for  a  particular  drag  coefficient  versus  Mach  table,  , 

,  and  s,.  are  unique  (independent  of  particle  mass,  air  density,  etc.).  Thus,  they  can  be 

precalculated  and  stored  in  tabulated  form. 


4.3  Mach  Number  as  a  Function  of  Scaled  Distance 

Equation  17  can  be  generalized  to  a  piecewise  function  of  scaled  distance: 

A 


M  = 


for  k  such  that  <s  <  . 


(36) 


Figure  7  presents  an  overview  of  the  functions  necessary  to  calculate  Mach  as  a  function  of 
scaled  distance.  Note  that  for  a  particular  drag  coefficient  versus  Mach  table,  a,^,  c^,  and 

s^.  are  unique  (independent  of  particle  mass,  air  density,  etc.).  Thus,  they  can  be  precalculated 

and  stored  in  tabulated  form. 
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4.4  Scaled  Time  as  a  Function  of  Mach  Number 

Equation  24  can  be  generalized  to  a  piecewise  function  of  Mach  number: 


1 


a 


T  =  ■ 


^In 

k 


a. 


M 


+  di. 


(37) 


for  k  such  that  <M  <  . 


Values  for  d^.  can  be  found  by  generalizing  equation  25: 


‘^k  ~  ^0,/t 


1  , 
- +  — ^In 

Pi 


f 


^k  + 


Pk 


M 


(38) 


0,k  J 


Values  for  initial  conditions  for  each  segment  of  the  function  are  given  by  equations  32  and  39: 


“  '^k 


(39) 


Values  for  can  be  obtained  from  equation  37: 


1 


a 


^k  = 


k-\ 


Pk-\^k  Pk-\ 


In 


Pk-. 


M 


+  d 


k-\ 


(40) 


k  J 


Since  the  value  of  is  arbitrary,  zero  is  used  as  a  convenient  choice: 


^-t=o  =  0 


(41) 


Figure  8  presents  an  overview  of  the  functions  necessary  to  calculate  scaled  time  as  a  function  of 
Mach  number.  Note  that  given  the  choice  of  =  0 ,  for  a  particular  drag  versus  Mach  table, 
a^,  p^,  d^,  and  are  unique.  Thus,  they  can  be  precalculated  and  stored  in  tabulated  form. 
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T 


+  dj^  for  k  such  that  <  M  < 


Figure  8.  Scaled  time,  r  ,  as  a  piecewise  function  of  Mach  number,  M. 


5.  State  Variables  and  Functional  Notation 


For  a  given  set  of  atmospheric  and  projectile  characteristics,  knowledge  of  any  one  of  the  three 
scaled  variables  defined  in  section  2  (i.e.,  M,  5,  or  r )  is  sufficient  to  uniquely  define  the  state  of 
a  projectile.  Applications  of  the  equations  presented  in  section  4  typically  involve  solving  for 
some  unknown  state  variable  given  one  or  two  known  state  variables.  The  following  six 
equations  can  make  solving  for  unknowns  easier. 

5.1  State  Variables 

Since  M,  s,  and  r  can  each  be  used  to  characterize  the  state  of  a  projectile,  the  following  three 
definitions  can  be  used  to  quantify  a  change  in  a  state  variable. 


initial 


(42) 


As  =  s 


’  _  o 

'final  ‘^initial 


(43) 


At  =  T 


—  7" 

final  initial 


(44) 


5.2  Functional  Notation 


Using  equation  30,  define  s{M)  to  be  scaled  distance  as  a  function  of  Mach  number: 


f,  M) 


(45) 


Using  equation  36,  define  M{s)  to  be  Mach  number  as  a  function  of  scaled  distance: 
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(46) 


M{s)  = 


A 


Using  equation  37,  define  r(M)  to  be  scaled  time  as  a  function  of  Mach  number: 

(  ^  ^ 

i  c 

V 


t(M)  =  — - -%ln 

A' 


+ 


A 

M 


+  ^/r 


(47) 


6.  Solving  for  Unknowns 


When  attempting  to  solve  for  some  unknown  state  variable,  it  helps  to  break  the  problem  into 
four  sections: 

•  Identify  known  and  unknown  values  for  unsealed  state  variables. 

•  Convert  known  values  for  state  variables  to  scaled  state  variables  using  the  equations 
presented  in  section  2. 

•  Use  the  equations  presented  in  section  5  to  develop  an  equation  that  is  solved  for  the 
unknown  scaled  state  variable  in  terms  of  known  scaled  state  variables. 

•  Use  the  equations  presented  in  section  2  to  convert  the  unknown  scaled  state  variable  to  an 
unsealed  value. 

For  example,  suppose  that  initial  speed,  ,  and  distance  traveled.  Ax ,  are  known,  and  it  is 
wished  to  determine  final  speed,  .  The  first  step  is  to  use  equations  1  and  2  to  convert 
and  Ax  to  the  scaled  variables  and  A^ .  Next,  represent  the  problem  in  functional 

notation: 

(48) 

where  M(.s^„^,)  is  given  by  equation  46. 

Solve  equation  42  for  the  final  scaled  speed  and  substitute  into  equation  46: 

M^„,,=M(A^  + 

^initial  ^  (49) 

Next,  note  that  can  be  calculated  using  equation  45.  Thus, 

+  (50) 

Finally,  use  equation  1  to  find  given  .  A  wide  variety  of  problems  can  be  solved  in 
this  manner. 
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7.  Example:  The  Lapua  Scenar  GB528  Rifle  Bullet 


7.1  Creating  a  Table  of  Unitless  Constants 

Wikipedia’s  External  Ballistics  page"^  lists  drag  coefficient  versus  Mach  number  values  for  the 
Lapua  Scenar  GB528  rifle  bullet.  Using  those  values,  equations  from  section  4  can  be  used  to 
calculate  the  parameters  listed  in  table  1 .  Specifically,  given  M and  ^  for  0  <  A:  <  28 ,  the 

equations  summarized  in  figures  7  and  8  have  been  used  to  calculate  ,  /?^ ,  5^ ,  ,  and 

df. . 


Table  1.  Unitless  parameters  for  a  Lapua  Scenar  GB528  rifle  bullet. 


k 

M 

C  D 

a 

p 

s 

c 

T 

d 

0 

2.4 

0.27 

-0.06 

0.414 

0 

5.2773 

0 

-0.241615 

1 

2.2 

0.282 

-0.05 

0.392 

0.315209 

5.555784 

0.137177 

-0.353937 

2 

2 

0.292 

-0.06 

0.412 

0.647242 

5.317506 

0.295484 

-0.237973 

3 

1.8 

0.304 

-0.085 

0.457 

1.000723 

4.89244 

0.481805 

-0.01001 

4 

1.6 

0.321 

-0.07 

0.433 

1.377521 

5.087262 

0.703845 

-0.139846 

5 

1.5 

0.328 

-0.08 

0.448 

1.576392 

4.96971 

0.832223 

-0.049923 

6 

1.4 

0.336 

-0.07 

0.434 

1.784183 

5.072469 

0.975621 

-0.13983 

7 

1.3 

0.343 

-0.05 

0.408 

2.002448 

5.268108 

1.137427 

-0.347739 

8 

1.2 

0.348 

0 

0.348 

2.234102 

5.791212 

1.32293 

-1.071706 

9 

1.15 

0.348 

0.02 

0.325 

2.3564 

6.034291 

1.427044 

-1.474873 

10 

1.1 

0.347 

0.08 

0.259 

2.48432 

6.938917 

1.54079 

-3.345151 

11 

1.075 

0.345 

0.16 

0.173 

2.550765 

9.120311 

1.601894 

-9.851053 

12 

1.05 

0.341 

0.28 

0.047 

2.619369 

26.548368 

1.666471 

-161.15269 

13 

1.025 

0.334 

1.12 

-0.814 

2.690775 

1.313248 

1.735305 

1.038474 

14 

1 

0.306 

2.8 

-2.494 

2.768002 

2.293195 

1.811601 

1.679499 

15 

0.975 

0.236 

2.36 

-2.065 

2.862002 

2.175026 

1.906853 

1.618416 

16 

0.95 

0.177 

0.92 

-0.697 

2.988737 

0.577959 

2.038621 

0.366765 

17 

0.925 

0.154 

0.48 

-0.29 

3.150184 

-3.032027 

2.210906 

-4.293853 

18 

0.9 

0.142 

0.2 

-0.038 

3.335448 

-45.258439 

2.413998 

-224.10354 

19 

0.875 

0.137 

0 

0.137 

3.537428 

17.072049 

2.641631 

-5.700392 

20 

0.85 

0.137 

-0.16 

0.273 

3.749016 

10.434934 

2.886984 

2.496055 

21 

0.825 

0.141 

-0.12 

0.24 

3.963785 

11.324716 

3.143445 

1.773405 

22 

0.8 

0.144 

-0.2 

0.304 

4.179723 

9.820507 

3.409242 

3.008442 

23 

0.7 

0.164 

-0.07 

0.213 

5.046777 

11.859988 

4.567078 

0.099253 

24 

0.6 

0.171 

-0.29 

0.345 

5.96672 

9.605173 

5.987225 

4.214716 

25 

0.5 

0.2 

-0.29 

0.345 

6.949257 

9.605173 

7.779309 

4.214716 

26 

0.4 

0.229 

-0.0025 

0.23 

7.988528 

10.413495 

10.102174 

-0.741033 

27 

0 

0.23 

- 

- 

- 

- 

- 

- 

All  of  the  parameters  listed  in  table  1,  though  specific  to  the  Lapua  Scenar  GB528  rifle  bullet, 
are  independent  of  air  density  and  initial  speed. 


^  Wikipedia.  External  Ballistics.  http://en.wikipedia.org/wiki/Extemal_ballistics  (accessed  3  June  2013). 
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7.2  Calculating  Speed  as  a  Function  of  Initial  Speed  and  Distance  Traveled 

For  the  Lapua  Scenar  GB528  rifle  bullet,  Wikipedia’s  External  Ballistics  page  lists  mass  {m)  as 
19.44  g,  diameter  as  8.59  mm,  and  initial  speed  (vinuiai)  as  830  m/s. 

Assuming  standard  sea-level  atmospheric  conditions  (air  density  (/?)  equals  1.225  kg/m^  and 
speed  of  sound  (v,)  equals  340.3  m/s)^  the  speed  of  a  bullet  that  has  traveled  300  m  ( Ax )  can  be 
calculated  using  the  equations  presented  in  sections  1  and  5. 


First  use  equation  1  to  find  initial  Mach  number  (Minuiai)- 


^initial  = 


^initial 
^  sound 


830  m/s 
340.3  m/s 


2.43902 


(51) 


Next,  use  the  equation  for  the  area  of  a  circle  to  find  the  cross-sectional  area  (A)  of  the  projectile: 


A  =  Tir^ 


3.14159 


^00859mY 
V  2  J 


=  5.7953x10"^  m" 


(52) 


Use  equation  2  to  find  the  scaled  distance  traveled  ( ): 


^  ^  ^  (1.225 kg/„»)(5.7953xl0-^n,^)  ^ 

2m  2(.01944kg) 


(53) 


Use  equation  45  to  find  the  initial  scaled  distance  (sinmai)'- 


1 


Pk=0 


In 


+ 


Pk=0 


M 


+  <^1=0  - 


1 


initial  J 


.414 


In 


-  .06  + 


V 


.414 
2.43902 


+  5.2773  =  -.059991  (54) 


Use  equation  43  to  find  the  final  scaled  distance  ): 


Sflnai  +  A?  = -.05999 1  +  .54778  =  .487789 


Use  equation  46  to  find  the  final  Mach  number  {M ): 

Pk^,  .392 


M 


c:  6c 


.392(.487789-5.555784) 


+  .05 


=  2.09454 


(55) 


(56) 


Finally,  use  equation  1  to  find  the  final  speed  ( 


Vflna,  =2.09454(340.3  m/s)  =  713  m/s 


(57) 


^  National  Aeronautics  and  Space  Administration.  U.S.  Standard  Atmosphere,  1976;  NASA-TM-X-74335;  U.S.  Government 
Printing  Office:  Washington,  DC,  October  1976. 


13 


7.3  Calculating  Time  of  Flight  as  a  Function  of  Initial  Speed  and  Distance  Traveled 

Use  equation  47  to  find  andr^„,, : 


1 


a 


k=0 


'  initial 


P, 


In 


A:=0 


P initial 

1 

(.414)(2.439) 

-.024768, 


a 


Pk=0 


k=0 


initial  J 


d 


yt=0 


(-.06) 


In 


-.06  + 


.414 

2.439 


+  (-.2416) 


final 


1 


a 


k=\ 


Pk=\l^  final 

1 


Pi 


In 


k=\ 


a 


k=\ 


Pk-l 

initial  J 


+  d 


k=\ 


^ _ (-.05) 

(.392)(2.095)  .392" 

=  .217245. 


In 


.05 


.392 

2.095 


+  (-.3539) 


Next,  calculate  the  change  in  scaled  time  ( Ar ): 

Ar  =  -  ^initial  =  -217245  -  (-.024768)  =  .242013 

Finally,  use  equation  3  to  find  time  of  flight  ( A/ ) 

2ot A T  _  2(.0 1 944  kg)(.2420 1 3) 


At  =  - 


pAv^  (1.225  kg/m  )(5. 795x10  m  )(340.3m/s) 


=  .390s 


(58) 


(59) 


(60) 


(61) 


8.  C++  Implementation 


8.1  Storing  Drag-Table  Coefficients:  The  DRAG  Struct 

DRAG  structs  are  used  to  store  the  parameters  associated  with  a  single  linear  drag  profile. 
Multiple  DRAG  structs  can  be  used  to  describe  a  nonlinear  drag  profile. 

DRAG  Code 


struct  DRAG{//<=============A  LINEAR  DRAG  PROFILE  (CREATE  A  SET  USING  Drags()) 

double  M;//< - THE  MAXIMUM  MACH  VALUE  FOR  THE  DRAG  PROFILE 

double  a;//< - THE  SLOPE  OF  THE  DRAG  PROFILE 

double  b;//< - THE  Y-INTERCEPT  OF  THE  DRAG  PROFILE 

double  Sj//< . . --THE  MINIMUM  SCALED  DISTANCE  FOR  THE  DRAG  PROFILE 

double  Cj//< . . . THE  SCALED-DISTANCE  CONSTANT  OF  INTEGRATION 

double  T;//< - THE  MINIMUM  SCALED  TIME  FOR  THE  DRAG  PROFILE 

double  d;//< - THE  SCALED-TIME  CONSTANT  OF  INTEGRATION 


}; //~~~YAGENAUT@GMAIL . COM  /N/ /N/ /N/ r\/ IN/ IN/ r\/ IN/ /N/ /N/ /N/ /N/ /N/ /N/ r\/ IN/ IN/ /N/ /N/ /N/ /N/ /N/ /N/ /N/ /N/  LAST-UPDATED-09DEC2013  IN/ IN/ IN/ IN/ IN/ IN/ 
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DRAG  Parameters 


M  M  specifies  the  maximum  Mach  value  for  the  drag  profile, 
a  a  specifies  the  slope,  a,  of  the  drag  profile  (see  equation  28). 

b  b  specifies  the  y-intercept,  p,  of  the  drag  profile  (see  equation  29). 

s  s  specifies  the  minimum  scaled  distance,  ,  for  the  drag  profile  (see  equations  33-35). 

c  c  specifies  c,  the  integration  constant  from  equation  3 1 . 

T  T  specifies  the  minimum  scaled  time,  r,  for  the  DRAG  profile  (see  equations  39^1). 

d  d  specifies  d,  the  integration  constant  from  equation  38. 

8.2  Creating  a  Pointer  to  an  Array  of  DRAG  Structs:  The  Drags()  Function 

The  Drags()  function  can  be  used  to  create  a  pointer  to  an  array  of  DRAG  structs.  The  code  uses 
equation  32  to  determine  ^ ;  equation  28  to  determine  ;  equation  29  to  determine  ; 

equations  33-35  to  determine  ^ ;  equation  31  to  determine  ;  equations  39-41  to  find  ; 
and  equation  38  to  determine  .  If  is  found  to  be  very  close  to  zero,  it  is  replaced  with  the 
user-definable  value  b_min. 

Note  that  the  Drags()  function  uses  the  “new”  command  to  allocate  memory  for  the  array  of 
DRAGS  that  is  pointed  to  by  the  return  value.  Thus,  to  avoid  memory  leaks,  each  use  of  the 
Drags()  function  should  be  accompanied  by  a  use  of  the  “delete[]”  operator. 

DragsQ  Code 

inline  DRAG*Drags(//<===================================CREATES  A  SET  OF  DRAGS 

const  double*M,//< - MACH  VALUESS  (MUST  BE  UNIQUE  &  DECREASING) 

const  double*Cd,//< . -CORRESPONDING  Cd  VALUES  (FROM  F=- . 5*rho*A*Cd*v''2) 

int  n,//< - THE  NUMBER  OF  MACH  VALUES 

double  b_min=lE-8){//< - MINIMUM  BETA  VALUE 

DRAG*D=new  DRAG[n-l]; 

fon(int  k=0j j=k-l j l«n-l;++k,++j){ 

D[k] .M=M[k]; 

D[k] .a=(Cd[k]-Cd[k+l])/(M[k]-M[k+l]); 

D[k] .b=Cd[k]-D[k] . a*M[k] ;/*<-*/if (fabs(D[k] . b)<b_min)D[k] .b=b_min; 

D[k] .s=!k?0:l/D[j] .b*log(D[j] .a+D[j] . b/M[k] )+D[ j ] .c; 

D[k] .c=D[k] .s-log(D[k] .a+D[k] . b/M[k] )/D[k] .b; 

D[k] .T=!k?0: (1/M[k]-D[j] .a/D[j] .b*log(D[j] .a+D[j] .b/M[k]))/D[j] .b+D[j] .d; 
D[k].d=D[k].T-(l/M[k]-D[k].a/D[k].b*log(D[k].a+D[k].b/M[k]))/D[k].b;} 

return  D;// . note  that  D  points  to  newly  allocated  memory 

}//—  -~YAGENAUT@GMAIL . COM  LAST-UPDATED-09DEC2013 
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DragsQ  Parameters 

M  M  points  to  an  array  of  Mach  values.  The  Mach  values  must  be  unique  and  listed  in 
decreasing  order. 

Cd  Cd  points  to  an  array  of  drag  coefficients.  There  must  be  a  drag  coefficient  for  each 
Mach  value. 

n  n  specifies  the  number  of  elements  contained  in  the  array  that  is  pointed  to  by  M. 

b_min  b_min  specifies  the  smallest  value  that  (3^  is  allowed  to  have  (see  the  discussion  that 
follows  equation  25). 

DragsQ  Return  Value 

The  DragsQ  function  returns  a  pointer  to  an  array  of  DRAG  structs. 

8.3  Calculating  Scaled  Distance  as  a  Function  of  Mach  Number:  The  SofMQ  Function 

The  SofMQ  function  uses  equation  45  to  calculate  scaled  distance  as  a  function  of  Mach  number. 

SofMQ  Code 


inline  double  SofM(//<============SCALED  DISTANCE  AS  A  FUNCTION  OF  MACH  NUMBER 

const  DRAG*D,//< - DRAG  PROFILES  (CREATE  USING  Drags()) 

int  n,//< - NUMBER  OF  MACH  VALUES  USED  TO  CREATE  D  (NOT  SIZE  OF  D!) 

double  M){//<--- - - - - - MACH  NUMBER  (M=v/c) 


int  k;/*<-*/for(k=n-2;k>0&&M>D[k] .M; --k); 

return  1/D[k]  .b*log(D[k]  .a+D[k]  .b/M)+D[k]  .c;// . s=rho*A/(2*m)*x 

}//—  -~YAGENAUT@GMAIL . COM  LAST-UPDATED-09DEC2013 

SofMQ  Parameters 

D  D  points  to  an  array  of  DRAG  structs.  The  CreateDragsQ  function  can  be  used  to 

create  a  pointer  to  an  array  of  DRAG  structs. 

n  n  specifies  the  number  of  Mach  elements  used  to  create  D.  Note  that  n  is  one  less 

than  the  size  of  D. 

M  M  specifies  a  Mach  number.  Use  equation  1  to  convert  from  speed  to  Mach  number. 

SofMQ  Return  Value 

The  SofMQ  function  returns  a  scaled  distance.  Use  equation  2  to  convert  from  scaled  distance  to 
measured  distance. 

8.4  Calculating  Mach  Number  as  a  Function  of  Scaled  Distance:  The  MofSQ  Function 

The  MofSQ  function  uses  equation  46  to  calculate  Mach  number  as  a  function  of  scaled  distance. 
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MofSO  Code 


inline  double  MofS(//<============MACH  NUMBER  AS  A  FUNCTION  OF  SCALED  DISTANCE 

const  DRAG*D,//<-- . . . ---DRAG  PROFILES  (CREATE  USING  Dnags()) 

int  n,//< - NUMBER  OF  MACH  VALUES  USED  TO  CREATE  D  (NOT  SIZE  OF  D!) 

double  s){//< - SCALED-DISTANCE  (s=nho*A/(2*m)*x) 


int  k; /*<-*/for(k=n-2; k>0&&s<D[k] . s; --k); 

return  D[k]  .b/(exp(D[k]  .b*(s-D[k].c))-D[k].a);// . M=v/c 

}//—  -~YAGENAUT@GMAIL . COM  /N/ /N/ /N/ r\/ r\/ IN/ IN/ /N/ /N/ /N/ /N/ /N/ /N/ r\/ IN/ IN/ /N/ /N/ /N/ /N/ /N/ /N/ /N/ /N^  LAST-UPDATED-09DEC2013  IN/ IN/ IN/ IN/ IN/ IN/ 


MofS()  Parameters 

D  D  points  to  an  array  of  DRAG  structs.  The  CreateDrags()  function  can  be  used  to 

create  a  pointer  to  an  array  of  DRAG  structs. 

n  n  specifies  the  number  of  Mach  elements  used  to  create  D.  Note  that  n  is  one  less 

than  the  size  of  D. 

s  s  specifies  a  scaled-distance.  Use  equation  2  to  convert  from  measured  distance  to 

scaled  distance. 

MofS()  Return  Value 

The  MofS()  function  returns  a  Mach  number.  Use  equation  1  to  convert  from  Mach  number  to 

speed. 

8.5  Calculating  Scaled  Time:  The  TofM()  Function 

The  TofM()  function  uses  equation  47  to  calculate  scaled  time  as  a  function  of  Mach  number. 

TofMO  Code 


inline  double  TofM(//<================SCALED  TIME  AS  A  FUNCTION  OF  MACH  NUMBER 

const  DRAG*D,//< - DRAG  PROFILES  (CREATE  USING  Drags()) 

int  n,//< - NUMBER  OF  MACH  VALUES  USED  TO  CREATE  D  (NOT  SIZE  OF  D!) 

double  M){//< - MACH  NUMBER  (M=v/c) 


int  k;/*<-*/for(k=n-2;k>0&&M>D[k] .M; --k);//  T=rho*A/(2*m)*c*t  - . 

return (1/M-D[k] .a/D[k] . b*log(D[k] . a+D[k] . b/M) )/D[k] . b+D[k] .dj//< . . ' 

}// - YAGENAUT@GMAIL  ■  COM - LAST~UPDATED~09DEC2013 - 


TofMQ  Parameters 

D  D  points  to  an  array  of  DRAG  structs.  The  CreateDrags()  function  can  be  used  to 

create  a  pointer  to  an  array  of  DRAG  structs. 

n  n  specifies  the  number  of  Mach  elements  used  to  create  D.  Note  that  n  is  one  less 

than  the  size  of  D. 

M  M  specifies  a  Mach  number.  Use  equation  1  to  convert  from  speed  to  Mach  number. 
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TofMQ  Return  Value 


The  TofM()  function  returns  a  scaled  time.  Use  equation  3  to  convert  from  scaled  time  to 
measured  time. 

8.6  Example:  Calculating  Speed  and  Time  of  Flight  for  the  Lapua  Scenar  GB528  Rifle 
Bullet 

The  following  example  uses  information  from  Wikipedia’s  External  Ballistics  page"^  to  calculate 
speed  and  time  of  flight  as  functions  of  total  distance  traveled  for  a  Lapua  Scenar  GB528  rifle 
bullet.  Calculated  values  are  compared  to  Doppler  radar  measurements,  which  were  also 
obtained  from  Wikipedia’s  External  Ballistics  page. 

ftinclude  <cstdio>// . printf() 

#include  "y_tna j_ld . h"// . yTna jlD, <cmath>{fabs( ) } 

int  main(){ 

double  M[]={2. 4, 2. 2, 2, 1.8, 1.6, 1.5, 1.4, 1.3, 1.2, 1.15, 1.1, 1.075, 1.050, 1.025,1, 

.975, .950, .925, .9, .875, .85, .825, .8, .7, .6, .5, .4,0}; 
double  C[]={. 270, .282, .292, .304, .321, .328, .336, .343, .348, .348, .347, .345, .341, 
.334, .306, .236, .177, .154, .142, . 137, . 137, . 141, . 144, . 164, . 171, . 2, . 229, . 23); 
yTnajlD: :DRAG*D=yTrajlD: :Dnags(M,C,28); 

double  A=3.14159*pow(8. 59/2000, 2);//. .cross-sectional  area  of  projectile  (m''2) 

double  m=19. 44/1000; // . mass  of  projectile  (kg) 

double  v0=830;// . initial  speed  of  projectile  (m/s) 

double  rho=1.225;// . density  of  air  at  sea  level  (kg/m''3) 

double  c=340.3;// . local  speed  of  sound  (m/s) 

double  v_r [ ]={830, 711, 604, 507,422, 349, 311,288, 267, 247, 227}; 
double  t_r[]={0, .3918, .8507,1.3937,2.0435,2.8276,3.748,4.7522,5.8254,7.0095, 
8.2909}; 

printf("%13s  |  calc.  ref.  |  calc.  ref.  \n",""); 

printfC'%13s  range  j  speed  speed  %%  |  time  time  %%\n",""); 

printf("%13s  (m)  |  (m/s)  ,  (m/s)  ,  diff.  |  (s)  ,  (s)  ,  diff .\n", ""); 

printf("%13s . \n",""); 

for(int  i=0; i<ll;++i){ 

double  v=yTrajlD: :MofS(D,28,yTrajlD: :SofM(D,28,v0/c)+300*i*rho*A/(2*m))*c; 
double  t=(yTrajlD: :TofM(D,28,v/c)-yTrajlD: :TofM(D,28,v0/c))*2*m/(rho*A*c); 
printf("%19.1f  |%6.1f  ,%6.1f  ,%5.1f  |%7.4f  ,%7.4f  ,%5.1f\n", 

300. *i, v, v_r[i] ,fabs(v_r[i] -v)/( (v_r [i]+v)/2)*100,t,t_r[i] , 
i==0?0:fabs(t_r[i] -t)/( (t_r[i]+t)/2)*100); } 
delete[]  D; 

}// - YAGENAUT@GMAIL .  COM - LAST~UPDATED~09DEC2013 - 
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OUTPUT: 


1 

calc . 

ref. 

1  calc. 

ref. 

range 

1 

speed 

speed 

% 

1  time 

time 

% 

(m) 

1 

(m/s) 

y 

(m/s) 

y 

diff . 

1  (s) 

y 

(S) 

y 

diff. 

0.0 

1 

830.0 

y 

830.0 

y 

0.0 

1  0 . 0000 

y 

0 . 0000 

y 

0.0 

300.0 

1 

712.7 

y 

711.0 

y 

0.2 

1  0.3901 

y 

0.3918 

y 

0.4 

600 . 0 

1 

606.3 

y 

604.0 

y 

0.4 

1  0.8467 

y 

0.8507 

y 

0.5 

900.0 

1 

509.3 

y 

507.0 

y 

0.5 

1  1.3867 

y 

1.3937 

y 

0.5 

1200.0 

1 

423.3 

y 

422.0 

y 

0.3 

1  2.0333 

y 

2.0435 

y 

0.5 

1500.0 

1 

350.2 

y 

349.0 

y 

0.3 

1  2.8140 

y 

2.8276 

y 

0.5 

1800.0 

1 

311.1 

y 

311.0 

y 

0.0 

1  3.7319 

y 

3.7480 

y 

0.4 

2100.0 

1 

288.3 

y 

288.0 

y 

0.1 

1  4.7344 

y 

4.7522 

y 

0.4 

2400 . 0 

1 

267.7 

y 

267.0 

y 

0.3 

1  5.8183 

y 

5.8254 

y 

0.1 

2700.0 

1 

245.2 

y 

247.0 

y 

0.7 

1  6.9894 

y 

7.0095 

y 

0.3 

3000.0 

1 

223.3 

y 

227. e 

y 

1.7 

1  8.2681 

y 

8.2909 

y 

0.3 

8.7  Example:  Calculating  Scaled  Time  of  Impact  for  Fragment  Trajectories 

The  following  example  simulates  the  trajectories  of  fragments  produced  by  a  fragmentation 
weapon. 

The  code  generates  1,000,000  fragments,  all  with  the  same  starting  location.  Each  fragment  is 
given  a  randomly  chosen  velocity.  Direction  is  uniformly  distributed  across  all  possible 
directions.  Speed  is  uniformly  distributed  with  values  from  Mach  1.0  to  Mach  3.0. 

A  rectangular  surface  that  is  80.0  x  40.0  scaled  length  units  is  placed  with  its  center  10.0  scaled 
length  units  away  from  the  fragment  source.  The  surface  is  oriented  such  that,  at  its  center,  it  is 
perpendicular  to  a  ray  originating  from  the  fragment  source. 

Distance  to  point  of  impact  is  calculated  using  functions  from  the  y3DOps  namespace.^  The 
MofSQ  function  is  used  to  find  each  fragment’s  speed  at  the  impact  location.  The  TofMQ 
function  is  used  to  find  the  scaled  time  of  impact. 


#include  <cstdio>// . f close () , FILE, f reopen () , printf ( ) , stdout 

#include  <cstdlib>// . rand( ) ,  RAND_MAX 

#include  "y_3d_ops . h"// . y3D0ps 

#include  "y_traj_ld.h"// . yTrajlD,<cmath>{cos()jSin() jsqrt()} 

int  main(){ 

FILE*f=freopen("impacts.txt", "w", stdout);// . redirect  output  to  a  file 

double  M[]={3,2.2,1.5,1,  .5,  .2},C[]={.45,  .5,  .6,  .9,  .5,  .2};// . Cd(M) 

yTrajlD: :DRAG*D=yTrajlD: :Drags(M,Cj6); 

double  T[]={10,40,20  ,  10,-40,20  ,  10,40,-20};// . rectangular  target 

forfint  i=0; i<1000000;++i){ 

double  M0=2*rand( )/RAND_MAX+l; // . initial  speed 

double  phi=2*3.14159265358979*rand()/RAND_MAX;// . random  azimuthal  angle 

double  theta=acos(2.*rand()/RAND_MAX-l);// . random  inclination  angle 


°Yager,  R.  J.  Three-Dimensional  Translations,  Rotations,  and  Intersections  Using  C++;  ARL-TN-557;  U.S.  Army  Research 
Laboratory:  Aberdeen  Proving  Ground,  MD,  2013. 
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double  S[ ]={0^0^0^ sin (theta )*cos( phi ) j sin (theta )*sin( phi ) ^ cos (theta ) }; 

double  t [3] ;/*<-*/y3D0ps : :IParameters(tjT^S);// . intersection  parameters 

double  x[3];/*<-*/y3DOps: :Intersect(XjtjS);// . point  of  intersection 

double  M=yTrajlD: :MofS(D,6,yTrajlD: :SofM(D,6,M0) 
+sqrt(x[0]*x[0]+x[l]*x[l]+x[2]*x[2] ) ) ; 

double  T=yTrajlD: :TofM(D^6jM)-yTrajlD: :TofM(D^6^M0);//. scaled  time  of  impact 
if(0<t[0]&&0<t[l]&&t[l]<l&&0<t[2]&&t[2]<l) 

printfC'%8.3f,%8.3f,%8.3f,%8.3f,%8.3f\n",x[0],x[l],x[2],M,T);} 
fclose(f ) ; 
delete[]  D; 

}//  /N/ /N/ /N/ /N/ /N/ IN/  YAGENAUT@GMAIL.COM  /N/ /N/ /N/ IN/ IN/ IN/ IN/ IN/ /N/ /N/ /N/ /N/ /N/ /N/ IN/ IN/ /N/ /N/ /N/ /N/ /N/ /N/ /N/ LAST-UPDATED-09DEC2013  IN/ IN/ IN/ IN/ IN/ IN/ 


The  results  are  shown  in  figure  9.  Since  all  time  and  length  units  are  scaled,  the  image  presented 
in  figure  9  is  independent  of  fragment  mass,  fragment  cross-sectional  area,  and  air  density.  Note 
the  color  mixing  that  results  from  selecting  random  initial  speeds. 


Figure  9.  Scaled  time  of  impact,  r  ,  for  fragments  striking  a  plate. 


9.  Code  Summary 


A  summary  sheet  is  provided  at  the  end  of  this  report.  It  presents  the  yTraj  ID  namespace,  which 
contains  the  struct  and  four  functions  that  are  described  in  sections  8. 1-8.5.  Also  presented  are 
the  examples  found  in  sections  8.6  and  8.7. 
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yTrajlD  Summary 


y traj ld.h 


#ifndef  Y_TRAD_1D_GUARD//  See  Yager,  R.j  Flanders,  B.  "A  One-Dimensional,  Nonit- 
#define  Y_TRAD_1D_GUARD//  erative  Trajectory  Model  (with  a  C++  Implementation)" 

#include  <cmath>// . exp(),fabs(),log() 

namespace  yTrajlD{//@(S)@@@@@@@@@@@@@@@^  A  ID,  NON-ITERATIVE  TRAIECTORY  MODEL 
struct  DRAG{//<=============A  LINEAR  DRAG  PROFILE  (CREATE  A  SET  USING  Drags()) 

double  Mj//< - THE  MAXIMUM  MACH  VALUE  FOR  THE  DRAG  PROFILE 

double  aj//< - THE  SLOPE  OF  THE  DRAG  PROFILE 

double  bj//< - THE  Y-INTERCEPT  OF  THE  DRAG  PROFILE 

double  Sj//< - THE  MINIMUM  SCALED  DISTANCE  FOR  THE  DRAG  PROFILE 

double  cj//< - THE  SCALED-DISTANCE  CONSTANT  OF  INTEGRATION 

double  Tj//< - THE  MINIMUM  SCALED  TIME  FOR  THE  DRAG  PROFILE 

double  dj//< - THE  SCALED-TIME  CONSTANT  OF  INTEGRATION 

}j// - YAGENAUT@GMAIL.COM - LAST~UPDATED~09DEC2013 - 

inline  DRAG*Drags(//<===================================CREATES  A  SET  OF  DRAGS 

const  double*M,//< - MACH  VALUESS  (MUST  BE  UNIQUE  &  DECREASING) 

const  double*Cd,//< - CORRESPONDING  Cd  VALUES  (FROM  F=- . 5*rho*A*Cd*v''2) 

int  n,//< - THE  NUMBER  OF  MACH  VALUES 

double  b_min=lE-8){//< - MINIMUM  BETA  VALUE 

DRAG*D=new  DRAG[n-l]j 

for (int  k=0, j=k-lj l«n-lj++k,++j){ 

D[k] .M=M[k]j 

D[k].a=(Cd[k]-Cd[k+l])/(M[k]-M[k+l])j 

D[k].b=Cd[k]-D[k].a*M[k]j/*<-*/if(fabs(D[k].b)<b_min)D[k].b=b_minj 

D[k].s=!k?0:l/D[j].b*log(D[j].a+D[j].b/M[k])+D[j].cj 

D[k].c=D[k].s-log(D[k].a+D[k].b/M[k])/D[k].bj 

D[k].T=!k?0:(l/M[k]-D[j].a/D[j].b*log(D[j].a+D[j].b/M[k]))/D[j].b+D[j].dj 

D[k].d=D[k].T-(l/M[k]-D[k].a/D[k].b*log(D[k].a+D[k].b/M[k]))/D[k].bj} 

return  Dj// . note  that  D  points  to  newly  allocated  memory 

}// - YAGENAUT@GMAIL.COM - LAST~UPDATED~09DEC2013 - 

inline  double  SofM(//<============SCALED  DISTANCE  AS  A  FUNCTION  OF  MACH  NUMBER 

const  DRAG*D,//< - DRAG  PROFILES  (CREATE  USING  Drags()) 

int  n,//< - NUMBER  OF  MACH  VALUES  USED  TO  CREATE  D  (NOT  SIZE  OF  D! ) 

double  M){//< - MACH  NUMBER  (M=v/c) 

int  kj/*<-*/for(k=n-2jk>0&&M>D[k] .Mj --k) j 

return  1/D[k]  .b*log(D[k]  .a+D[k]  .b/M)+D[k]  .cj// . s=rho*A/(2*m)*x 

}// - YAGENAUT@GMAIL.COM - LAST~UPDATED~09DEC2013 - 

inline  double  MofS(//<============MACH  NUMBER  AS  A  FUNCTION  OF  SCALED  DISTANCE 

const  DRAG*D,//< - DRAG  PROFILES  (CREATE  USING  Drags()) 

int  n,//< - NUMBER  OF  MACH  VALUES  USED  TO  CREATE  D  (NOT  SIZE  OF  D! ) 

double  s){//< - SCALED-DISTANCE  (s=rho*A/(2*m)*x) 

int  kj/*<-*/for(k=n-2jk>0&&s<D[k] .sj --k) j 

return  D[k]  .b/(exp(D[k]  .b*(s-D[k]  .c))-D[k]  .a)j// . M=v/c 

}// - YAGENAUT@GMAIL.COM - - - LAST~UPDATED~09DEC2013 - 

inline  double  TofM(//<================SCALED  TIME  AS  A  FUNCTION  OF  MACH  NUMBER 

const  DRAG*D,//< - DRAG  PROFILES  (CREATE  USING  Drags()) 

int  n,//< - NUMBER  OF  MACH  VALUES  USED  TO  CREATE  D  (NOT  SIZE  OF  D! ) 

double  M){//< - MACH  NUMBER  (M=v/c) 

int  kj/*<-*/for(k=n-2jk>0&&M>D[k] .Mj --k) j//  T=rho*A/(2*m) *c*t  - . 

return(l/M-D[k] .a/D[k] .b*log(D[k] .a+D[k] .b/M))/D[k] .b+D[k] .dj//< - ' 

}// - YAGENAUT@GMAIL.COM - LAST~UPDATED~09DEC2013 - 


"Impacts"  Example 


#include  <cstdio>// . fclose(),  FILE,f  reopen()  ,printf  (),  stdout 

#include  <cstdlib>// . rand( ) ,  RAND_MAX 

#include  "y_3d_ops.h"// . y3D0ps 

#include  "y_traj_ld.h"// . yTrajlD,<cmath>{cos(),sin(),sqrt()} 

int  main(){ 

FILE*f=freopen("impacts.txt", "w",stdout)j// . redirect  output  to  a  file 

double  M[]={3,2.2,1.5,1,  .5,  .2},C[]={.45,  .5,  .6,  .9,  .5,  .2}j// . Cd(M) 

yTrajlD: :DRAG*D=yTrajlD: :Drags(M,C,6) j 

double  T[]={10,40,20  ,  10,-40,20  ,  10,40, -20}j// . rectangular  target 

for(int  i=0ji<1000000j++i){ 

double  M0=2*rand()/RAND_MAX+lj// . initial  speed 

double  phi=2*3.14159265358979*rand()/RAND_MAXj// . random  azimuthal  angle 

double  theta=acos(2.*rand()/RAND_MAX-l)j// . random  inclination  angle 

double  S[ ] ={0,0,0, sin (theta) *cos(phi), sin(theta)*sin (phi), cos(theta)}j 

double  t[3] j/*<-*/y3D0ps: : IParameters(t, T, S) j // . intersection  parameters 

double  x[3] j/*<-*/y3D0ps: :Intersect(x,t,S)j// . point  of  intersection 

double  M=yTrajlD: :MofS(D,6,yTrajlD: :SofM(D,6,M0) 
+sqrt(x[0]*x[0]+x[l]*x[l]+x[2]*x[2] )) J 

double  T=yTrajlD: :TofM(D,6,M)-yTrajlD: :TofM(D,6,M0)j//.scaled  time  of  impact 
if(0<t[0]&&0<t[l]&&t[l]<l&&0<t[2]&&t[2]<l) 

printfr'%8.3f,%8.3f,%8.3f,%8.3f,%8.3f\n",x[0],x[l],x[2],M,T)j} 
fclose(f) J 
delete[]  Dj 

}// - YAGENAUT@GMAIL.COM~— — - - - LAST~UPDATED~09DEC2013 - 


Graph  of  Output  from  "Impacts"  Example 


Equations 


"Validating"  Example 


#include  <cstdio>// . printf() 

#include  "y_tra j_ld .h"// . yTrajlD,<cmath>{fabs() } 

int  main(){ 

double  M[] ={2.4, 2. 2, 2, 1.8, 1.6, 1.5, 1.4, 1.3, 1.2, 1.15, 1.1, 1.075, 1.050, 1.025,1, 
.975, .950, .925, .9, .875, .85, .825, .8, .7, .6, .5, .4,0}j 
double  C[]={.270, .282, .292, . 304, . 321, . 328, .336, .343, .348, .348, .347, .345, .341, 
.334, .306, .236, .177, .154, . 142, . 137, . 137, .141, .144, .164, .171, .2, .229, .23}j 
yTrajlD: :DRAG*D=yTrajlD: :Drags(M,C,28) j 

double  A=3.14159*pow(8.59/2000,2)j//. .cross-sectional  area  of  projectile  (m''2) 

double  m=19.44/1000j// . mass  of  projectile  (kg) 

double  v0=830j// . initial  speed  of  projectile  (m/s) 

double  rho=1.225j// . density  of  air  at  sea  level  (kg/m"3) 

double  c=340.3j// . local  speed  of  sound  (m/s) 

double  v_r [ ]={830,711,604, 507,422,349, 311, 288, 267, 247, 227} j 
double  t_r[]={0, .3918, .8507,1.3937,2.0435,2.8276,3.748,4.7522,5.8254,7.0095, 
8.2909}j 

printf("%13s  |  calc.  ref.  |  calc.  ref.  \n","")j 

printf("%13s  range  j  speed  speed  %%  |  time  time  %%\n","")j 

printf("%13s  (m)  j  (m/s)  ,  (m/s)  ,  diff.  |  (s)  ,  (s)  ,  diff . \n", "" ) j 

printf ( "%13s======================================================\n " , " " ) j 

for(int  i=0ji<llj++i){ 

double  v=yTrajlD: :MofS(D, 28, yTrajlD: :SofM(D, 28,v0/c)+300*i*rho*A/(2*m) )*cj 
double  t=(yTrajlD: :TofM(D,28,v/c)-yTrajlD: :TofM(D, 28,v0/c) )*2*m/(rho*A*c) j 
printf ("%19. If  |%6.1f  ,%6.1f  ,%5.1f  |%7.4f  ,%7.4f  ,%5.1f\n", 

300. *i, v,v_r [i],fabs(v_r [i] -v)/((v_r[i]+v)/2)*100,t,t_r[i], 
i==0?0:fabs(t_r[i] -t)/((t_r[i]+t)/2)*100) j } 
delete[]  Dj 

}// - YAGENAUT@GMAIL.COM - - - - - LAST~UPDATED~09DEC2013 - 


Output  from  "Validating"  Example 


range 

(m) 

calc, 
speed 
(m/s)  , 

ref. 
speed 
,  (m/s)  , 

% 

diff. 

1  calc. 

1  time 

1  (s)  . 

ref. 

time 

,  (s)  , 

% 

,  diff. 

0.0 

830.0  , 

,  830.0  , 

0.0 

1  0.0000  , 

,  0.0000  , 

,  0.0 

300.0 

712.8  J 

,  711.0  , 

0.2 

1  0.3901  , 

,  0.3918  , 

,  0.4 

600.0 

606.1  J 

,  604.0  , 

0.3 

1  0.8466  , 

,  0.8507  , 

,  0.5 

900.0 

509.3  J 

,  507.0  , 

0.5 

j  1.3867  , 

,  1.3937  , 

,  0.5 

1200.0 

423.2  , 

,  422.0  , 

0.3 

1  2.0333  , 

,  2.0435  , 

,  0.5 

1500.0 

350.2  J 

,  349.0  , 

0.3 

j  2.8140  , 

,  2.8276  , 

,  0.5 

1800.0 

311.2  J 

,  311.0  , 

0.0 

1  3.7319  , 

,  3.7480  , 

,  0.4 

2100.0 

288.2  , 

,  288.0  , 

0.1 

1  4.7344  , 

,  4.7522  , 

,  0.4 

2400 . 0 

266.6  J 

,  267.0  , 

0.1 

j  5.8163  , 

,  5.8254  , 

,  0.2 

2700.0 

245.1  J 

,  247.0  , 

0.8 

1  6.9894  , 

,  7.0095  , 

,  0.3 

3000.0 

224.0  J 

,  227.0  , 

1.3 

j  8.2697  , 

,  8.2909  , 

,  0.3 
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